Hofstadter's butterfly in moire superlattices: A fractal quantum HaU effect 



C. R. Deani'2, L. Wang^, P. Maher^ C. Forsythe^ F. Ghahari^ Y. Gao^, J. Katoch"^, M. Ishigami"^, 
P. Moon^ M. Koshino^ T. Taniguchi^, K. Watanabe^, K. L. Shepard\ J. Hone^ and P. Kim^ 

^Department of Electrical Engineering, Columbia University, New York, NY 
^Department of Mechanical Engineering, Columbia University, New York, NY 
^Department of Physics, Columbia University, New York, NY 
"^Department of Physics and Nanoscience Technology Center, University of Central Florida, Orlando, FL 
^ Department of Physics, Tohoku University, Sendai, Japan and 
^National Institute for Materials Science, 1-1 Namiki, Tsukuba, Japan 



PACS numbers: xx 

Electrons moving through a spatially periodic 
lattice potential develop a quantized energy spec- 
trum consisting of discrete Bloch bands. In two 
dimensions, electrons moving through a magnetic 
field also develop a quantized energy spectrum, 
consisting of highly degenerate Landau energy 
levels. In 1976 Douglas Hofstadter theoretically 
considered the intersection of these two prob- 
lems and discovered that 2D electrons subjected 
to both a magnetic field and a periodic electro- 
static potential exhibit a self-similar recursive en- 
ergy spectrum^. Known as Hofstadter's butter- 
fly, this complex spectrum results from a delicate 
interplay between the characteristic lengths as- 
sociated with the two quantizing fields^ and 
represents one of the first quantum fractals dis- 
covered in physics. In the decades since, exper- 
imental attempts to study this effect have been 
limited by difficulties in reconciling the two length 
scales. Typical crystalline systems (< 1 nm peri- 
odicity) require impossibly large magnetic fields 
to reach the commensurability condition, while 
in artificially engineered structures (> 100 nm), 
the corresponding fields are too small to com- 
pletely overcome disorde r^^^^^l Here we demon- 
strate that moire superlattices arising in bilayer 
graphene coupled to hexagonal boron nitride pro- 
vide a nearly ideal-sized periodic modulation, en- 
abling unprecedented experimental access to the 
fractal spectrum. We confirm that quantum Hall 
effect features associated with the fractal gaps 
are described by two integer topological quantum 
numbers, and report evidence of their recursive 
structure. Observation of Hofstadter's spectrum 
in graphene provides the further opportunity to 
investigate emergent behaviour within a fractal 
energy landscape in a system with tunable inter- 
nal degrees of freedom. 

The total number of electron states per area of a com- 
pletely filled Bloch band is exactly no = where A 
is the area of the unit cell of the periodic potential. In 
a magnetic field, the number of states per area of each 
filled Landau level (LL) is given hy B/(j)Q where = h/e 
is the magnetic fiux quanta. The quantum description 



of electrons subjected simultaneously to both a periodic 
electric field and a magnetic field can be simply param- 
eterized by the dimensionless ratio where (j) = BA 
is magnetic fiux per unit cell. The general solution of 
this problem, however, exhibits a rich complexity due 
to the incommensurate periodicities between the Bloch 
and Landau states^^. In his seminal work^, Hofstadter 
showed that for commensurate fields, corresponding to 
rational values of = p/q^ where p and q are co- 
prime integers, the single-particle Bloch band splits into 
q subbands (beginning with the Landau level descrip- 
tion it is equivalently shown^ that at these same rational 
values each Landau level splits into p subbands). This 
results in a quasi-continuous distribution of incommen- 
surate quantum states that exhibits self-similar recursive 
structure, yielding the butterfiy-like fractal energy dia- 
gram (see SI). 

Important insight into this problem was subsequently 
provided by Wanniei^, who considered the density of 
charge carriers, n, required to fill each fractal subband. 
Replotting the Hofstadter energy spectrum as integrated 
density versus field, Wannier realized that all spectral 
gaps are constrained to linear trajectories in the density- 
field diagram. This can be described by a simple Dio- 
phantine relation 

{n/uo) = t{(l)/(l)o) + s (1) 

where n/uo and (^/(po are the normalized carrier density 
and magnetic fiux, respectively, and s and t are both 
integer valued. Here n/rio represents the Bloch band 
filling fraction, which is distinct from the usual Landau 
level filling fraction, u = tk^q/B (the two are related by 
the normalized fiux, i.e. n/rio = i^^/^o)- The physi- 
cal significance of the s and t quantum numbers became 
fully apparent with the discovery of the integer quantum 
Hah effect^^ in 1980, after which it was shown, simulta- 
neously by both Stfreda^ and Thouless et al^, that the 
Hall conductivity associated with each mini-gap in the 
fractal spectrum is quantized according to (Jxy = te'^/h. 
The second quantum number, s, physically corresponds 
to the Bloch band filling index in the fractal spectrunP. 
This formalism suggests several unique and unambiguous 
experimental signatures associated with the Hofstadter 
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energy spectrum that are distinct from the conventional 
quantum Hall effect: i) the Hall conductance can vary 
non-monotonically and can even fluctuate in sign, ii) Hall 
conductance plateaus together with vanishing longitudi- 
nal resistance can appear at non-integer LL filling frac- 
tions, iii) the Hall conductance plateau remains quan- 
tized to integral multiples of e^//i, however, the quanti- 
zation integer is not directly associated with the usual LL 
filling fraction. Instead, quantization is equal to the slope 
of the gap trajectory in the n/n o versus (j) / (j)o Wannier 
diagram, in accordance with the Diophantine equation. 

Mini-gaps within the fractal energy spectrum becom e 
significant only once the magnetic length {Ib = ^Jh/eB)^ 
which characterizes the cyclotron motion, is of the same 
order as the wavelength of the periodic potential, which 
characterizes the Bloch waves. For usual crystal lattices, 
where the inter-atomic spacing is a few angstroms, the 
necessary magnetic field is impossibly large, in excess 
of 10,000 T. The main experimental effort therefore has 
been to lithographically define artificial superlatticeJP''^ 
with unit cell dimension of order tens of nanometers so 
that the critical magnetic field remains small enough to 
be achievable in the lab, yet still large enough so that 
the quantum Hall effect is fully resolved without being 
smeared out by disorder. Fabricating the optimally-sized 
periodic lattice, while maintaining coherent registry over 
the full device and without introducing substantial dis- 
order has proven a formidable technical challenge. Pat- 
terned GaAs/AlGaAs heterostructures with ~ 100 nm 
periodic gates provided the first experimental support 
for Hofstadter's predict ionJ^^^nH. However, limited abil- 
ity to tune the carrier density or reach the fully devel- 
oped quantum Hall effect (QHE) regime in these sam- 
ples has made it difficult to map out the complete spec- 
trum. While similar concepts have also been pursued in 
non-solid-state model systems^^ the rich physics of the 
Hofstadter spectrum remains largely unexplored. 

Heterostructures consisting of atomically thin materi- 
als in a multi-layer stack provide a new route towards 
realizing a two-dimensional system with laterally mod- 
ulated periodic structure. In particular, coupling be- 
tween graphene and hexagonal boron nitride (hBN), 
whose crystal lattices are isomorphic, results in a periodic 
moire pattern. The moire wavelength is directly related 
to the angular rotation between the two lattice^^^*^, 
and is tunable through the desired length scales with- 
out the need for lithographic techniques^ Moreover 
hBN provides an ideal substrate for achieving high mo- 
bility graphene devices, crucial for high resolution quan- 
tum Hall measurements^^ while field effect gating in 
graphene allows the Fermi energy to be continuously var- 
ied through the entire moire Bloch band. 

In this study, Bernal-stacked bilayer graphene (BLG) 
Hall bars are fabricated on hBN substrates (Fig. [l^,b) 
using mechanical exfoliation followed by co-lamination 
(see methods) . Fig. lb shows a non-contact atomic force 
microscopy (AFM) image acquired from an example de- 
vice. In the magnified region, a moire pattern is visible 
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FIG. 1. Moire super lattice (a) Cartoon schematic of 
graphene on hBN showing the emergence of a moire pat- 
tern. The moire wavelength varies with the mismatch angle, 
0. (b) Left shows an AFM image of a multi-terminal Hall 
bar. Right shows a high resolution image in a magnified re- 
gion. The moire pattern is evident as a triangular lattice (up- 
per inset shows a further magnified region). FFT of the scan 
area (lower inset) confirms a triangular lattice symmetry with 
period 15.5 =b 0.9 nm. (c) Resistance measured versus gate 
voltage at zero magnetic field. Inset shows the corresponding 
conductivity versus temperature, indicating that the satellite 
features disappear above ^100 K. (d) Longitudinal resistance 
(left axis) and Hall resistance (right axis) versus gate volt- 
age dX B — 1 T. The Hall resistance inverts sign and passes 
through zero at the same gate voltage as the satellite peaks. 



with triangular symmetry. Fast Fourier transform (FFT) 
analysis of the image, shown in the inset, indicates that 
the moire wavelength is 15.5 ± 0.9 nm. This is compa- 
rable to the maximal wavelength of ~ 14 nm expected 
for graphene on hB]>l^^^^ (set by the 1.8% lattice mis- 
match between the two crystals), suggesting that in this 
device the BLG lattice is oriented with near zero angle 
mismatch to the underlying hBN lattice. 

Fig. ^ shows transport data measured from the same 
device. In addition to the usual resistance peak at the 
charge neutrality point (CNP), occurring at gate voltage 
Vg 2 \^ two additional satellite resistance peaks ap- 
pear, symmetrically located at Vsati ~ ±30 V away from 
the CNP. These satellite features are consistent with a 
depression in the density of states (DOS) at the super- 
lattice Brillouin zone band edge, analogous to previous 
spectroscopic measurements of single layer graphene cou- 
pled to a moire pattern^^^^. Assuming non-overlapping 
bands, |Vsatz|, gives an estimate of the moire wavelength 
to be ~ 14.6 nm (see supplemental information), in good 
agreement with the AFM measurements. The nature of 
these satellite peaks can be further probed in the semi- 
classical, low 5-field transport regime. In Fig. [T]l, Ion- 
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FIG. 2. Emergence of anomalous quantum Hall states (a) Landau fan diagrams showing longitudinal resistance, Rxx, 
(left) and Hall resistance, Rxy (right). Rxx is plotted against magnetic field on the vertical and gate bias on the horizontal. 
In Rxy the axes are scaled by the size of the moire unit cell to give ^/^o on the vertical and n/rio on the horizontal. QHE 
states corresponding to the conventional BLG spectrum are indicated by white lines. Solid yellow and red lines track QHE 
outside this conventional spectrum, with dashed lines indicating the projected n/noo intercept. The slope of each line is labeled 
across the top axis, (b) Longitudinal and transverse Hall conductivities corresponding to line cuts at constant magnetic field 
(constant 0/00 ) from the Landau fan diagram in (a). At B = 7 T the QHE ladder is consistent with previous reports for 
bilayer graphene. At B = 18 T and 26 T additional QHE states emerge showing Hall conductivity plateaus quantized to integer 
multiples of e^//i, but appearing at non-integer LL filling fractions. Yellow and Red colored bars indicate correspondence to 
the same anomalous features marked by solid lines in (a). Blue bars label the conventional QHE features. Numbers label the 
quantization integer for each plateau. 



gitudinal resistance, R^x^ and transverse Hall resistance, 
Rxy, are plotted versus gate voltage at 5 = 1 T. Near 
the central CNP, the Hall resistance changes sign as the 
Fermi energy passes from the electron to the hole band. 
The same trend also appears near Vsath consistent with 
the Fermi energy passing through a second band edge. 
This provides further confirmation that the moire pat- 
tern, acting as a periodic potential superlattice, gives rise 
to a mini-Brillouin zone band^^. The satellite peaks van- 
ish at temperatures above 100 K (inset of Fig. Ic), 
indicating that the coupling between the BLG and hBN 
atomic lattices is of order ~10 meV. 

In the remainder of this letter, we focus on magneto- 
transport measured at high field. Fig. [2^ shows the evo- 
lution of Rxx and Rxy for magnetic fields up to 31 T. In 
the left panel Rxx is plotted against the experimentally 
tunable gate voltage and magnetic field (a so-called Lan- 
dau fan diagram). In the right panel, the magnitude of 
the corresponding Rxy is plotted against the dimension- 
less parameters appearing in the Diophantine equation, 
n/uo and (/)/0o- This "Wannier diagram" is simply the 
Landau fan diagram with both axes relabeled by dimen- 
sionless units defined by normalizing to the moire unit 
cell area. 

In a conventional quantum Hall system, the Landau 
fan diagram exhibits straight lines, tracking minima in 
Rxx and plateaus in Rxy Plotted against n/uo and (l)/(l)oi 



the slope of each line is precisely the LL filling fraction, 
and all lines converge to the origin. White lines in 
Fig. [2^ identify QHE states matching this description, 
tracking LL filling fractions u = 4^8 and 12. This is 
consistent with the usual QHE hierarchy associated with 
a conventional degenerate BLG spectrum. 

At large magnetic fields, several additional QHE states, 
exhibiting minima in Rxx together with plateaus in Rxy, 
develop outside the usual BLG sequence that also follow 
straight lines in the Landau fan diagram but converge to 
non-zero values of n/uo. Yellow and red lines in Fig. [2|i 
trace examples of these anomalous QHE states appear- 
ing within the lowest LL. Unlike the conventional QHE 
states, each of the anomalous QHE states is character- 
ized by both an integer valued intercept, s (yellow and 
red lines converge to n/n^ = 1 and 2, respectively) and 
an integer valued slope, t (labeled along the top axis in 
the figure). In Fig. [2]3, longitudinal and Hall conductivi- 
ties acquired at constant magnetic field (corresponding to 
horizontal line cuts through the fan diagram in Fig. |2|l) 
are plotted against LL filling fraction, u. At large mag- 
netic fields the anomalous QHE states are remarkably 
well developed, exhibiting wide plateaus in (Jxy concomi- 
tant with zero valued axx- Moreover, these states appear 
in general at non-integer filling fractions. Comparison be- 
tween Fig. [2^ and (2)3 further reveals that Hall conductiv- 
ity plateaus are quantized to integer values te^ /h^ where 
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FIG. 3. Fractal gaps (a) Landau fan diagrams similar to Fig. 2 but measured from a separate device. Here the zero- field 
satellite peak position indicates a moire period of 11.6 nm or approximately 1.5 x smaller superlattice unit cell area. Significantly 
more structure is observed compared with the previous device (b). Lower panel shows the evolution of Gxx with temperature 
varying between 2 K and 20 K, acquired at constant B — 2^ T corresponding to the line cut shown in (a). Upper panel shows 
the corresponding Gxy at T = 2 K. The bracketed numbers label the (s,t) values of the corresponding fractal gaps according 
to the Diophantine equation, (c) Bubble plot of energy gaps determined from the temperature dependence acquired at two 
magnetic fields = 25 T and 28.5 T). The gaps are plotted as circles with radius scaled relative to the largest gap value 
measured. Dashed lines trace out select fractal gap positions allowed by the Diophantine equation. Solid lines trace regions 
where the corresponding fractal gaps appear as minimia in Gxx together with quantized plateaus in Gxy 



the quantization integer t equals the slope in the Wannier 
diagram. Similar internal structure is observed within 
higher order Landau levels (see Fig. [s] and also supple- 
mental info). The anomalous QHE states observed here 
are consistent with fully developed spectral gaps result- 
ing from a Hofstadter-type energy spectrum. Moreover, 
our ability to fully map the density-field space provides 
a remarkable confirmation of the Diophantine equation, 
where we observe direct evidence that QHE features asso- 
ciated with the Hofstadter spectral gaps are characterized 
by the two quantum numbers, s, and t, corresponding to 
the n/rio intercept and slope, respectively, in the Wannier 
diagram. 

Fig. |3] shows similar data to Fig. |2j but acquired 
from a separate device in which the moire wavelength 
is only 11.6 nm. Again, QHE states appear outside of 
the conventional Bernal-stacked BLG sequence following 
straight lines whose origin and slope are both integer val- 
ued, and with the slope exactly matching the Hall quanti- 
zation, in precise agreement with the Diophantine equa- 
tion. Similar to the previous device, the v = {) insulating 
state undergoes a dramatic change near (/)/(/)o = 1/2 when 
anomalous QHE states associated with the fractal gaps 
begin to fully develop. 

In Fig. [sJd the lower panel shows cjxx acquired at 
5 = 25 T, corresponding to a horizontal line cut through 
Fig. [3^, with varying sample temperature. The magni- 
tude of the fractal gaps were estimated from the tempera- 
ture dependence of the cr^x minima in the thermally acti- 
vated regime (see methods). This analysis was performed 
at two magnetic field values, 5 = 25 T and 28.5 T. Fig.jsj) 
summarizes our findings. Each fractal gap is marked by 



a circle centered at the corresponding n/no, (j)/(j)o co- 
ordinate, and with radius scaled relative to the largest 
gap value (all gap values are given in the supplemen- 
tal information). With increasing magnetic field spectral 
energies develop in a complicated way where some gaps 
grow with field {e.g. (s^t) = (1,1) and (2,-2)), while 
others diminish {e.g. (1, —4)). At fixed magnetic field it 
appears generally true that for constant s values, fractal 
gap states with increasing t exhibit increasing gap size. 
For example, at 5 = 25 T, A(i _4) ~ 48 K, whereas 
A(i _3) ^ 30 K. This contradicts Hofstadter 's predic- 
tion that fractal gaps corresponding to lower quantum 
numbers exhibit larger gap values. We note that such a 
trend was subsequently found to be specific to square lat- 
tice symmetries^ ^ ^. Furthermore, a non-trivial case also 
arises when two fractal gap states overlajP, such as in our 
data between the (1,-2) and (2,-3) states as 0/0o ^ 1- 
Further theoretical analysis specific to moire-patterned 
BLG is necessary to fully understand the trends high- 
lighted here. 

Fig. [4^ shows a normalized Landau fan diagram of Rxy 
values corresponding to the same Rxx data in Fig. [3^. 
Dashed horizontal lines in the figure label special values 
(/)/(/)o = 1/m, where m is integer valued. Referred to by 
Hofstadter as the "pure cases" ^ , these lines of high sym- 
metry provide the framework for the recursive structure 
of the butterfly spectrum, marking the boundaries of the 
repeating sub-cells that appear within the main cell^i It 
appears that at particular points in the fan diagram, Rxy 
tends towards zero, inverting sign as the magnetic field 
is swept through these pure fractions. This is further re- 
vealed by the single line-trace in in Fig. ^p. At low mag- 
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FIG. 4. Recursive structure (a) Rxy Wannier diagram for the same device shown in Fig. 3. White sohd Unes label (j)/(j)o 
values corresponding to "pure fractions", (j)/(j)o — (b) Longitudinal conductivity, Gxx (left axis) and Hall conductivity, Gxy 

(right axis) versus magnetic field, acquired at a constant gate voltage {Vg = —39 V, corresponding to the white dashed line in 
(a). Blue, green and purple bands mark the boundaries of the conventional A/" = 2,3 and 4 Landau levels, respectively, where 
N is the index of the orbital wavefunction. Dashed lines mark the B field corresponding to the pure fractions labeled in (a), (c) 
Expanded plot showing details of both Gxx and Gxy in the vicinity of (/)/(/)o = 1/2 from (b). Plateaus in Gxy concomitant with 
mmima m Cxx resemble a a mini QHE trace centered around (t)/(j)o — 1/2, consistent with Hofstadter's prediction of recursion 
in the butterfly spectrum. Inset shows the position of the mini QHE plateaus plotted versus effective magnetic field following 
a linear trend (see text) 



netic fields, we observe the semi-classical Weiss oscilla- 
tions previously reported for square superlattices^^ 
In the quantum Hall regime, the longitudinal conduc- 
tivity shows a local peak as the magnetic field passes 
through the pure fractions, with the corresponding Hall 
conductivity exhibiting a sharp transition. Near the 
field corresponding to (j)/(l)o = 1/2, labeled in Fig. ^ 
as ^1/2, plateaus appear in R^y together with minima 
in Rxx resembling a mini QHE series centered around 
Bi/2- If we redefine the local effective magnetic field as 
B' = 5 — ^1/2, then according to the usual QHE for- 
malism we expect the relation = (1/ B')n'h/ where 
v' is an effective filling fraction given by the Hall quan- 
tization, B' labels the effective magnetic field position of 
the Rxx minima, and n' represents an effective carrier 
density. Inset in Fig. [Ij^ shows a plot of v' versus 1/B' 
and indeed, the data follows a linear trend. In spite of 
the large magnetic field (5i/2 ~ 17.3 T), this remarkable 
observation indicates that locally the electrons behave 



as if the magnetic field is reduced to zero. We regard 
this as compelling evidence of the long-predicted recur- 
sive nature of the Hofstadter spectrum where repeated 
mini fan diagrams emerge within the main one. Interest- 
ingly, the linear trend shown inset in Fig.|4j3 does not pass 
through the origin, but is vertically offset by 4.1 ± 0.1. 
The origin of this offset is unclear but may be related to 
disorder effects since in this regime the spectrum is not 
fully gapped^. 

Finally, we discuss the apparent broken symmetries ob- 
served in our data. We first note that the most prominent 
fractal gap states within the lowest Landau level of both 
devices correspond to positive 5, but negative t quantum 
numbers. Likewise, the overall fractal gap structure ap- 
pears stronger on the hole side (negative density values) 
than the electron side (positive density values). This elec- 
tron/hole asymmetry results from coupling between BLG 
and hBN, which, in a single-particle tight binding calcu- 
lation, tends to break the bipartite nature of the BLG 
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lattice. We further note that our data exhibits complete 
symmetry breaking of both spin and pseudospin degen- 
eracy, with fractal gaps appearing for both even and odd 
values of s and t, whereas only multiples of 4 may be 
expected in a fully degenerate case. The pseudospin de- 
generacy, which is related to the spacial inversion symme- 
try, is lifted due to asymmetric coupling to the substrate, 
where the bottom graphene layer interacts more strongly 
with the BN substrate than the top graphene layer. In 
the high magnetic field regime, Zeeman coupling may 
break spin symmetry to yield odd integer QHE states. 
A detailed consideration of these effects can be found in 
the supplemental information. However, we note that 
at these fields the Coulomb energy is more than an or- 
der of magnitude larger than Zeeman^^, suggesting that 
interaction-induced spontaneous symmetry breaking of 
the internal quantum degrees of freedom may be neces- 
sary to understand our experimental observations. 



METHODS 

Devices were fabricated using a co-lamination me- 
chanical transfer technique similar to that described 
previousl^ ^— I Electrical leads consisting of a 
Cr/Pd/Au metal stack were deposited using standard 
electron-beam lithography after which the sample was 
etched into a Hall bar by exposure to oxygen plasma. 
Graphene/hBN stacks were fabricated on doped Si sub- 
strates with a ~ 300 nm oxide layer. More than 20 de- 
vices were made in similar way, where 6 devices show 
similar behavior to that reported here. We focus only on 
2 high quality devices in the text, listing other examples 
in the supplementary materials. In electrical measure- 



ments the charge carrier density was varied using the 
doped silicon as a field effect gate. Four- terminal trans- 
port measurements were performed using a lock-in am- 
plifier at 17 Hz with a 10 — 100 nA source current. Sam- 
ples were measured in a 31 T resistive magnet and ^He 
cryostat (sample in vapour). Longitudinal and Hall con- 
ductivities were calculated from the corresponding mea- 
sured resistance according to = Pxx/{pxx ^xy) 
(Txy = RxyKplx ^Rly)^ respectively. 

AFM images were acquired after device fabrication was 
complete, using an Omicron low temperature AFM sys- 
tem operated at room temperature. Imaging was per- 
formed using Vi)ias = 0.2 V and 5f = 20 Hz. Images 
were filtered to remove noise. 

Gap energies were estimated from the temperature 
dependence of longitudinal conductivity minima in the 
thermally activated regime, a^x ex e^/^^^-^ where A is 
the energy gap, /c^ is Boltzmann's constant and T is the 
electron temperature. Each gap value was determined 
from the corresponding Arrhenius plot by fitting to the 
linear regime. 
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SUPPLEMENTARY INFORMATION 



SI. WANNIER DIAGRAM 



The relevant quantities in the Landau fan diagram for 
a 2D electron system with a superimposed periodic po- 
tential are the size of the superlattice unit cell, and the 
number of flux penetrating this unit cell. For a square 
superlattice the magnetic flux penetrating the unit cell is 





Bo? 



FIG. SI. Left is a cartoon image of an electron subjected 
to both a magnetic field, and a square periodic lattice. The 
usual Hofstadter Butterfly spectrum, calculated for a square 
superlattice, is shown in the middle. Right shows the spec- 
trum replotted as density versus magnetic field. Coloured 
lines in the Hofstadter energy-field spectrum (middle) follow 
lines of constant chemical potential. In the Wannier density- 
field spectrum (right), these same lines follow linear trajecto- 
ries, according to the Diophantine equation. 



Fig. [Sip and [Sip illustrates the relationship between 
the energy-field diagram and density-field diagram. As 
first demonstrated by Wannier (see main text), all spec- 
tral gaps follow linear trajectories in the density- field 
space according to a dimensionless Diophantine equation 
(Eqn. 1 in the main text). The Wannier diagram is ex- 
perimentally accessible by performing transport measure- 
ment while varying the carrier density and magnetic field. 
In the quantum Hall regime, spectral gaps in the spec- 
trum appear as minima in longitudinal resistance, Rxx^ 
and quantized plateaus in the transverse Hall resistance. 



S2. REPLOTTING LANDAU FAN DIAGRAM IN 
DIMENSIONLESS UNITS 



(2) 



where, B is the applied magnetic field, o? is the unit cell 
area, and (j)o = h/ e \^ the magnetic flux quantum. For a 

hexagonal superlattice, the unit cell area is ^o? where 
a is the unit cell basis vector, giving 



2h 



(3) 



We also normalize the electron density by the inverse area 
of the unit cell, i.e., n n/uo^ where 



1 



unit cell area 



(4) 



We can define the critical field as the field were (j)/<j)o ^ 
1, i.e. when the magnetic length approaches the scale 
were there is approximately 1 flux quantum per unit cell 
area. From [2] and [Sj this gives the corresponding critical 
magnetic field 



Bo 



2h 



(5) 



The same analysis for a square lattice gives Bq = h/ (ea^). 
Fig. [S2] shows the critical field versus unit cell wave- 
length. The ideal unit cell length is highlighted in the 
grey shaded region, spanning the magnetic field range 
from a few Tesla up to approximately 40 Tesla. 



S3. MOIRE SUPERLATTICES IN BILAYER 
GRAPHENE ON HBN 



hexagonal lattice 
square lattice 




20 40 60 80 100 120 140 160 
unit cell wavelength (nm) 



FIG. S2. Critical field versus unit cell length, defined where 
the ratio 0/0o 1. Gray shaded area highlights the ideal 
unit cell length for magnetic fields achievable in the labora- 
tory. 



Fig. S3 shows three different bilayer graphene/hBN 
devices that exhibit modulated transport due to the pres- 
ence of a moire pattern. Device si (Fig.[S3^) corresponds 
to the device in Fig. 2 of the main text, with the longest 
moire wavelength observed in our devices, s3 (Fig. [SS]^) 
corresponds to Fig. 3 of the main text, with the short- 
est moire wavelength amongst our devices. Device s2 
(Fig. S3 3) shows a device with an intermediate moire 
wavelength. Assuming a filled band model the moire 
wavelength was determined from the satellite peak po- 
sition according to 



rio 



9s9v 



(6) 
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FIG. S3, (a)-(c) shows transport measured at zero magnetic 
field (top panel) and corresponding Rxx plotted in a Landau 
fan diagram (bottom panel), for three different devices, which 
we label here S1,S2 and s3. In each panel, a indicates the cor- 
responding moire wavelength, determined from the position 
of the satellite resistance peaks at zero field (see text). 



where Usat is the field-effect density at the satellite peak 
position, rio = l/A where A = \fZ(j? 12 is the unit cell 
area of the moire pattern (assuming hexagonal symme- 
try; see Fig. 1 in the main text) with a the moire wave- 
length, is the electron spin degeneracy and is the 
valley, or pseudospin, degeneracy. Usat is calculated ac- 
cording to a standard parallel plate capacitor model to 
be Usat = Cgiy^at - VcNp)/e, where Cg is the geometric 
capacitance, Vsat is the gate voltage corresponding to the 
position of the satellite resistance peak and Vcnp is the 
gate voltage corresponding the central charge neutrality 
peak. Substituting back into equation [6] and solving for 
the moire wavelength gives 



V?>Cg{Vsat - Vcnp) 



(7) 



For the three devices in Fig|S3[ the satellite peaks ap- 
pear at 30 V, ,37.5 V, and 48.5 V, giving moire wave- 
lengths 14.6 nm, 13.1 nm and 11.6 nm, respectively. 



S4. THEORETICAL SPECTRUM FOR BILAYER 
GRAPHENE ON HBN 

The Diophantine equation represents a universal set of 
constraints describing all possible gaps within the Hofs- 
tadter spectrum. Which gaps remain open for a particu- 
lar device depends critically on the details of the system. 
Fig. S4a shows the calculated Hofstadter spectrum for 
Bernal-stacked bilayer graphene oriented to hBN. Here 
graphene and hBN are modeled by tight-binding honey- 
comb lattices with the lattice periods a ~ 0.246 nm and 
<^hBN ^ 0.2504 nm,^ respectively. We assume that Bernal 
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FIG. S4. (a) Energy spectrum calculated for bernal-stacked 
bilayer graphene on hBN with a 0° mismatch angle. The 
energy scale is chosen to highlight the lowest Landau level. 
The numbers label the quantized Hall conductivity in units 
/h corresponding to each gap. (b) Experimental magne- 
toresistance measured for device with similar angular mis- 
match (data is the same as shown in Fig. 2 in the main text). 
Yellow and red lines correspond to s = 1 and s = 2, respec- 
tively. Numbers along top axis label corresponding t values, 
(c) Wannier diagram calculated for the energy spectrum in 

(a) , showing only gaps within the lowest Landau level. Hori- 
zontal and vertical axes are density and magnetic flux, respec- 
tively, normalized to the area of the moire unit cell. Fractal 
gaps are plotted as black circles where the radius is propor- 
tional to the gap width. Coloured lines highlight gaps that are 
observed in the experimental data, (d) Similar calculation to 

(b) but including Zeeman coupling and interlayer asymmetry 
(see text) 



graphene bilayer and hBN monolayer are aligned with 
zero rotation angle, and the ratio between the two lattice 
constants is round to a rational number ahBN/<^ = 56/55 
to give a finite moire superlattice period 56a ^ 13.8 nm. 
We assume the interlayer distance of bilayer grpahene 
to be 0.335 nm and that between hBN and graphene 
0.322 nm^. We consider only a p orbital state on each 
atomic site, and set on-site potential to 0, 3.34 eV and 
— 1.40 eV for C, B and N atoms, respectively.l^ For the 
hopping amplitudes between different sites, we adopt the 
Slater-Koster parametrization used for twisted bilayer 
graphen^, irrespective of combination of atoms. To com- 
pute the low-energy spectrum in magnetic field, we take 
the low- lying Landau levels (|£^| < 1.5 eV) of isolated bi- 
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layer graphene as the basis P and the couphng with hBN 
states in the high-energy region is included as on-site po- 
tential on the bottom graphene layer (faced to hBN) in 
the second-order perturbation. 

In Fig. S4b, the energy spectrum is replotted as a 
Wannier type diagram where the energy axis is replaced 
by the normalized density. Fractal gaps inside the low- 
est Landau levels (between v = ±4) are plotted as black 
circles with radius proportional to the gap size. The the- 
oretical calculation correctly predicts the experimentally 
observed asymmetry in the lowest Landau level (see also 
Fig. 3 in the main text) where fractal gaps originat- 
ing from s > are much stronger than s < 0. Re- 
lated, the calculation correctly identifies a partial lift- 
ing of the 4- fold degeneracy in graphene, evident as fea- 
tures associated with s = ±2, rather than only multiples 
of 5 = 4. This results from broken valley degeneracy 
due to inversion-asymmetric coupling where the bottom 
graphene layer interacts more strongly with the BN sub- 
strate than the top graphene layer. Several details in 
the theoretical spectrum disagree with experimental ob- 
servations. For example, strong mini-gaps, predicted on 
the electron side, are not apparent in the experiment. 
Additionally, in experiment all symmetries are broken at 
high field, with fractal gaps corresponding to both odd 
and even integer values of s and t appearing, whereas the 
theoretical spectrum remains two-fold degenerate at all 
fields since the spin Zeeman splitting is neglected. 

In Fig. S4c, the Wannier diagram is recalculated in- 
cluding spin Zeeman splitting {g = 2) and interlayer po- 
tential difference between the top and the bottom layer 
of bilayer graphene. The latter is assumed to be the sum 
of the constant term 30 meV and the gate-induced term 
proportional to electron density with the relative dielec- 
tric constant Sr = To get a better agreement, we also 
magnified coupling strength between graphene and hBN 
by a factor about 1.4 compared to the original (Fig. S4c). 
While additional symmetry breaking states appear, ex- 
act correspondence with our data is difficult to achieve. 
The model depends sensitively on physical parameters, 
such as interlayer coupling strength and interlayer po- 
tential difference, that are not well known. Moreover, 



this model does not account for many-body interaction, 
which may play a substantial role at these large magnetic 
fields. Further theoretical and experimental work is nec- 
essary to better understand the details of the observerd 
fractal energy spectrum. 

S5. TEMPERATURE DEPENDENCE 



(a) B =25T 
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FIG. S5. Temperature dependent transport acquired at (a) 
B = 25 T and (b) B = 28.5 T. In both (a) and (b) left panel 
shows the longitudinal conductivity with temperature varying 
between 2 K and 20 K, together with Hall conductivity at 
2 K. Right panel shows Arrhenius plots for each of the fractal 
gaps labeled in the left panel. The corresponding gap value, 
determined from the slope in the linear regime, is given in 
each plot. 
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